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We present simulations of flux-gradient-driven superconducting vortices interacting with strong 
columnar pinning defects as an external field H(t) is quasi-statically swept from zero through a 
matching field B<f,. We analyze several measurable quantities, including the local flux density 
B(x,y,H(t)), magnetization M(H(t)), critical current J c (B(t)), and the individual vortex flow 
paths. We find a significant change in the behavior of these quantities as the local flux density 
crosses B^,, and quantify it for many microscopic pinning parameters. Further, we find that for a 
given pin density J C (B) can be enhanced by maximizing the distance between the pins for B < B^. 

PACS numbers: 74.60.Ge, 74.60.Jg 



Introduction. — Flux pinning by correlated disorder in 
superconductors has recently been the subject of intense 
investigations because of its ability to immobilize vor- 
tices, reduce dissipation effects, and create high critical 
currents Q. In particular, the introduction of correlated 
disorder by heavy ion irradiation can create columnar 
defects in which vortices are pinned along their entire 
length H. Such systems are predicted to form a Bose 
glass at low temperatures when the vortices localize at 
randomly-placed columnar pins The behavior of 

this state depends on whether the magnetic flux density 
B is below, equal to, or greater than B^, the matching 
field at which the number of flux lines equals the number 
of columnar pins In addition, an important role is 

played by the critical state |^||, in which a gradient in 
the flux profile arises as individual vortices enter the sam- 
ple and become pinned at defects. This gradient deter- 
mines J C (B) |l]JH|. Detailed information on the vortex 
dynamics in the critical state, including local and aver- 
aged physical quantities, is required to explain the effects 
induced by the addition of columnar defects. Further- 
more, understanding the effect of the spatial distribution 
of the pinning sites on the Bose glass can facilitate creat- 
ing samples with optimal pinning to enhance J c , which 
is of technological importance. To investigate the flux- 
gradient-driven dynamics of vortices in superconductors 
with columnar defects, we perform simulations of indi- 
vidual vortices and antivortices interacting with strong 
columnar pinning sites, as an external field H(t) is swept 
through a complete loop. We monitor distinct changes 
in relevant physical quantities (e.g., B, M, J c ). 

Simulation. — We simulate an infinite slab with a mag- 
netic field H = Hz applied parallel to the surface so that 
there are no demagnetization effects. We consider rigid 
vortices and straight columnar pins of uniform strength 
(all || z); thus we only need to model a transverse 2D 
slice (x-y plane) from the 3D slab. Our results are 
for a 36A x 36A system with periodic boundary condi- 



tions, where A is the penetration depth. Our super- 
conducting system has a pinned region from x = 6A 
to x — 30A; thus 1/3 of the system is unpinned and 
2 /3 has randomly- placed non-overlapping parabolic traps 
with radius £ p = 0.15A. The flux lines evolve according 
to a T = molecular dynamics algorithm. Thus, thermal 
effects are neglected and we consider a situation deep in 
the Bose glass regime, where the Mott insulator phase 
should be observable |ji||§. Note that pinning due to 
columnar disorder is much less sensitive to thermal ef- 
fects than pinning due to point defects. For instance, for 
B < B^ Tauber et al. [Q], working with Bi2Sr2CaCu20s, 
find that thermal effects on pinning become relevant only 
for T\ w 0.9T C w 78 K m Ti rrev . Therefore, in this case 
all temperatures T < T\ may be considered low. 

We simulate an increasing external field H(t) as de- 
scribed in |u|[llj where we quasi-statically add flux lines 
to the unpinned region; these flux lines attain a uniform 
density n, so that an external applied field H may be de- 
fined as H = n<I>o. As the external field increases, the flux 
lines — by their own mutual repulsion — are forced into the 
pinned region where their motion is impeded by the de- 
fects. We correctly model the vortex-vortex force by a 
modified Bessel function K\(r/X), which falls off expo- 
nentially so that a cutoff at 6A can be safely imposed. 
The overdamped equation of motion is: fi = f ™ + f" p = 
f]Vi, where the total force fi on vortex i (due to other vor- 
tices f™, and pinning sites f" p ) is fi = J^f=i fo ^i(l r i — 
r.'I/A^ + E^! t l r ^ r i P) l e((£ P -|r,-ri p) |)/A) r ik . 
Here O is the Heaviside step function, (vj) is the lo- 
cation (velocity) of the ith vortex, r£ p ^ is the location of 
the fcth pinning site, £ p is the radius of the pinning well, 
fp is the maximum pinning force, and we take t] = 1. 
We measure lengths in units of the penetration depth 
A, forces in terms of f — <&q/87t 2 A 3 (= 2e /A), ener- 
gies in units of eo = (<E>o/47rA) 2 (= /oA/2), and magnetic 
fields in units of <I>o/A 2 . The sign of the interaction be- 
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tween vortices is determined by /„; we take f v = +/o 
for repulsive vortex- vortex interactions and f v = —fo for 
attractive vortex-antivortex interactions. A vortex and 
antivortex annihilate and are removed from the system 
if they come within 0.3A of one another. For colum- 
nar pinning the condensation energy is eo = ($o/47rA) 2 , 
which can be related to the maximum pinning force f p 
to give f p / fo = A/£ p = A/0.15A « 7 . However the mea- 
surable quantities of the system do not change for pin- 
ning forces greater than 2.0/o- We systematically vary 
one parameter such as the density of pins, n p , pinning 
strength f p , or the spatial distribution of pinning, while 
keeping the other parameters fixed. Previous simulations 
of flux-gradient-driven vortices [|l^,[ll| considered situa- 
tions where B < and f p /fo < 1 and found a J C (H) 
dependence as 1/H. This corresponds to the experimen- 
tally common situation described by the Kim model 
However, when we consider very strong columnar pins 
(fp/fo > 2.0) we find a very different behavior. 

Flux Profiles. — From the location of the vortices, 
we compute their local density B(x,y, H(t)) and profile 
B(x,H(t)) inside the sample. Figure 1 shows a typical 
example of B{x) for a hysteresis loop. When \B\ < B^, 
the gradient of the flux profile (i.e., the current J c ) is 
very large. Above B^ we see an abrupt change in the 
flux profile to a much shallower slope. As the external 
field is ramped still higher, the slope of B decreases grad- 
ually. These results are very similar to those observed 
experimentally [|7||l3|] for samples with columnar pins as 
the local field is increased beyond B^ and quite differ- 
ent from the ones seen in other types of samples (e.g., 
|nj[nj). Such results have been interpreted in terms of 
a change in the pinning force: when B < B^ the vor- 
tices are pinned at individual defects, and when B > B^ 
many vortices are pinned at interstitial sites due to the 
repulsion from the defect-pinned vortices. 

To give a more microscopic idea of the behavior in this 
system, in Fig. 2(a-c) we present snapshots of the posi- 
tions of the pins and the vortices in an 8A x 8A region 
of the sample (while the external field is gradually in- 
creased). The vortices are entering from the left of the fig- 
ures and thus there is a flux gradient in the three frames. 
In Fig. 2(a), most of the vortices are pinned at defects; 
only two (located near other defect-pinned flux lines) are 
at interstitial sites. Since most of the vortices at these 
low fields are trapped by defects, the effective pinning 
force is very high and hence, as seen in Fig. 1, the flux 
gradient for low fields is quite steep. In Fig. 2(b), near 
B<p, a greater fraction of interstitially pinned vortices is 
present, producing a much lower flux gradient. The fact 
that unoccupied pins are typically located very close to 
occupied pins indicates that many unoccupied pins are 
being "screened" by the vortex repulsion from the oc- 
cupied pins. Also, because of the flux gradient, there 
are more interstitially pinned vortices near the left of the 
figure than the right. Thus, the presence of a flux gra- 



dient limits the number of accessible pinning sites, when 
these sites are randomly distributed. Finally in Fig. 2(c), 
above B^, the flux gradient is small and hence the effec- 
tive pinning force is low. Here a large number of vortices 
are weakly trapped at interstitial sites and the major- 
ity of the defects are now occupied by vortices. In the 
lower left corner a small domain of triangularly-arranged 
vortices can be seen — indicating that the vortex-vortex 
repulsion is starting to dominate over the vortex-pin at- 
traction. 

Magnetization loops and Critical Current. — From the 
microscopic vortex dynamics of individual vortices de- 
scribed above, we can compute macroscopic measurable 
quantities. In particular, we quantify the vortex behav- 
ior described above, as the field is brought through B^, 
by measuring M(H) and J C {B) for samples with differ- 
ent pinning parameters (see Fig. (3)). For instance, the 
effect of changing the matching field B^ is examined in 
Fig. 3(a,b) by fixing f p at 2.5/o and changing the pin 
density n p . For all the cases we considered, the width of 
the magnetization curve is much broader when \B\ < B^, 
falling off very rapidly when B > B$. Fig. 3 shows the 
enhancement of J C (B), for B < B^, as the field is swept. 
We obtain J C (B) directly from the flux density profiles 
using Maxwell's equation dB/dx — /j,qJ. For a fixed 
field H , we average the slope of B over the entire sam- 
ple. In all cases the enhancement of J C (B) is restricted 
to \B\ < B&. This is very similar to the results seen in 
p3[ at low T, where J C (B) was also obtained directly 
from the flux gradient. It can also be seen from Fig. 3 
that J c has a definite dependence on B, which is differ- 
ent from the step-like two-current model suggested in jTj 
(indicated by a dashed line in Fig. 3(b)). 

We next examine the effect of varying the pinning force 
from f p — 0.3/o to f p — 5.0/o (see Fig. 3(c,d)). There 
is little change in M(H) and J C (B) as the pinning force 
is increased from f p — 2.5/o to f p — 5.0/o- This reflects 
the fact that once vortices are trapped strongly enough, 
so that they cannot be unpinned, any further increase in 
the pinning strength does not significantly affect M(H) 
and J C (B). As the pinning force is decreased in Fig. 3(c), 
the width of the magnetization loop narrows. A crossover 
at Be/, is no longer observed when the pinning force from 
the individual defects becomes on the order of the inter- 
stitial pinning force (when f p ~ fo/2). This behavior is 
emphasized in Fig. 3(d) for J C (B). At the lowest pin- 
ning strengths, only a small increase in J C (B) can be 
seen when \B\ < B$/3, and not when \B\ < B^ as in 
the strongly pinned (f p > 2.0/o) cases. These results are 
very similar to those seen in E| when the temperature 
was raised, which reduces the effective pinning strength. 

The theory of a Mott insulator phase in the Bose glass 
H leads one to expect an enhancement in J C (B) at B^. 
At this field the vortices become strongly localized at 
the pinning sites, due to the repulsion of the other vor- 
tices, and produce a Meissner-like phase which prevents 
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increased flux penetration over a range of field. However, 
we do not observe any particular enhancement of J c at 
B = £>0. This is due to the combined effects of both 
the critical state, which prevents a uniform field in the 
sample, and the randomness in the spatial distribution 
of defects, which prevents all the pins from being occu- 
pied due to "screening" effects. We point out that near 
-B^ there can still be a decrease in the magnetization re- 
laxation rate since it depends on thermal activation and 
tunneling. 

We can significantly reduce screening effects (due to 
closely-spaced defects) by changing the spatial arrange- 
ment of the pinning sites, thus enhancing J c . Figure 
3(e,f) presents M{H) and J C (B) for three systems with 
pinning sites which are (1) randomly located, (2) placed 
in a triangular lattice (see, e.g., 10), and (3) randomly 
displaced, up to 0.25A from a triangular lattice. All three 
systems have B$ = 1/A 2 and f p — 2.5/q. For fields less 
then Be/,, the width of M(H) and J C {B) for the trian- 
gular pinning is nearly twice that for the random case. 
The asymmetry in M(H) and J C (B) is much more pro- 
nounced when pins are placed in a triangular lattice. For 
this case, there is also a very sharp drop in M(H) just 
past very similar to that observed in |L4}| , since here 
the distinction between vortex pinning at defects and 
interstitial pinning is much better defined. This result 
is very important since it indicates that substantial in- 
creases in the critical current for fields less then B^ can 
be achieved by evenly spacing the pins. In this case, 
M(H) and J C (B) are enhanced because all the pins are 
accessible to the vortices, as opposed to the random case 
where it is difficult for vortices to reach closely-spaced 
defects. Moreover, in the system with triangular pin- 
ning, trapped vortices are better localized at the defects 
due to the interactions with other pinned vortices (since 
an additional energy minima is present when the vor- 
tices are in a favorable lattice configuration). This indi- 
cates that the Mott insulator phase is more accessible in 
a triangular array of pinning sites than in a Bose glass 
with randommly-located defects. In the case where the 
pinning sites are randomly displaced (up to 0.25A) from 
the triangular sites, there is still a clear enhancement of 
J C (B); however, it is much smaller than in the ordered 
triangular pinning array. 

Vortex Plastic Flow. — To investigate the dynamics 
of the flux-gradient-driven vortices near B$ and higher 
fields, we present in Fig. 4(a,b) the trajectories of the 
vortices for a system with (a) strong pinning, f p — 2.5/o, 
and with (b) weaker pinning, f p — 0.3/oj both with 
£?0 = I.O^o/A 2 . For both panels the external field H 
is increased from 0.9<f>o/A 2 to 1.4$o/A 2 . In (a) the mov- 
ing vortices follow specific winding channels indicative of 
plastic flow of interstitial vortices around domains of flux 
lines which are strongly pinned at defects. These "vor- 
tex rivers" often develop well defined tight bottlenecks as 
well as much broader "vortex streets" with several lanes. 



At the top left side of Fig. 4(a), two bottlenecks can be 
seen. Each of these "interstitial vortex paths" is located 
between flux lines that are strongly pinned at defects, in- 
dicating that the interstitial vortices are flowing through 
the energy minima created by the strongly pinned flux 
lines. These interstitial-vortex-channels form for both 
the entrance and exit of vortices and are a general fea- 
ture of all simulations with randomly placed strong pin- 
ning sites. It is interesting to note that although a large 
number of vortices are clearly immobile at the defects, 
J c is quite low for this region of B > B^ indicating that 
the channels are acting as weak links. In Fig. 3(b) with 
f p = 0.3/o we see a much different behavior with only a 
few channels. Further, it can be observed from the lines 
that cross through pinning sites that vortices can become 
unpinned since the pinning force is on the order of the 
interstitial pinning. Thus vortex motion consists of both 
interstitial trajectories and also pin-to-pin motion. 

Conclusions. — Our simulations of flux-gradient-driven 
rigid vortices interacting with columnar pinning sites are 
consistent with recent experiments jjj which indicate that 
a sharp change in the magnetic flux gradient occurs above 
B$. Our results quantitatively predict how J c varies with 
field and pinning parameters (e.g., strength and loca- 
tion). Moreover, we can monitor the spatio-temporal dy- 
namics of vortices as the pinning mechanism evolves from 
strong columnar pinning (for B < B^) to weaker intersti- 
tial pinning (for B > B<p). In particular, we quantify the 
enhancement of J c , when B < B<p, by varying the spa- 
tial arrangements of the pinning sites. When the pinning 
strength of the columnar defects becomes on the order of 
the interstitial pinning strength, this sharp transition in 
J c is greatly reduced. We compute the spatio-temporal 
dynamics of vortices and show that for strongly pinned 
samples the vortex transport is dominated by the plas- 
tic flow of interstitially-pinned vortices around regions 
of flux lines strongly pinned at defects. In the weaker 
pinning samples, vortices can jump from defect to defect 
site. Finally, M and J c can be considerably enhanced by 
placing the defects on a regular triangular lattice, for a 
given pin density, so that all pins are equally accessible, 
thus preventing the screening of closely-spaced defects 
occurring in a random distribution. 
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FIG. 1. Magnetic flux density profiles, B(x,H(t)) = 

(36A) -1 J™ X dyB(x,y,H(t)), for a 36A x 36A system 
with a 24A x 36A region (the actual sample) contain- 
ing 864 pinning sites with / p =2.5/o, £ P =0.15A. Here 
B =1.O $ /A 2 =864$o/(36A x 24A). For the initial ramp-up 
phase in (1) a total of 2600 vortices are added, so the max- 
imum external field is 2.0 $ /A 2 ~2600$ /(36A x 36A). In 
(2) the field is ramped down and then reversed to a final 
value of — 2.0<&o/A 2 . Finally in (3), the field is brought up to 
2.0$o/A 2 . A large gradient in B can be seen for \B\ < B^. 
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FIG. 2. An 8A x 8A region of the 24A x 36A sample studied 
in Fig. 1. The vortices enter from the left edge of each frame 
(which corresponds to the sample edge). Panels (a), (b), and 
(c) correspond to B/B<p « 0.3, 1.0, and 1.6, respectively, of 
the ramp up phase in Fig. 1. Pinning sites are indicated by 
open circles, while vortices are shown as filled dots. 
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FIG. 3. Magnetization loops (top panels) and the corre- 
sponding critical currents (bottom panels) for several sam- 
ples. The J c (B)s are taken directly from the B(x) during 
"ramp-down" (e.g., stage (2) in Fig. 1). In (a,b) f p is held 
fixed at 2.5/o and the density of pinning sites n v is varied: 
B^/Qo = 0.50/A 2 , 0.75/A 2 , 1.0/A 2 . In (c,d) n v remains fixed 
(i?^ = 1.0<E>o/A 2 ) and the pinning strength f p is changed. In 
(e,f) B,t, = 1.0<I>o/A 2 , f p — 2.5/o, and the location of the pin- 
ning sites is varied, (f) shows the significant enhancement of 
J C (B) that results from defects placed in a regular triangular 
array, as opposed to random placement. These results show 
that even a distorted triangular array of pinning sites signifi- 
cantly enhances J C (B) over the case with a random location 
of pinning sites. 




(a) (b) 

FIG. 4. (a,b) show the trajectories of vortices while the 
external field is raised from 0.95 $o/A 2 to 1.4$o/A 2 . (a) A 
10A x 10A region of the sample used in Fig. 1. The strength 
of the pinning force, f p , is 2.5/o, for the strong-pinning case 
(a), and 0.3/o, for the weak-pinning case (b). In (a) the vor- 
tex transport is characterized by vortex trails of interstitial 
vortices which move around regions with flux lines that are 
strongly-pinned at defects. In (b) vortex transport proceeds 
in a different manner: pin-to-pin vortex motion, as well as in- 
terstitial, is possible and the previously-narrow vortex trails 
become considerably broader. 
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